13 Population-level models 
of forest dynamics 


M. D. Korzukhin and M. Ya. Antonovski 


13.1. Introduction 


It has been a long-standing tradition in mathematical ecology to use 
dynamical equations in the modeling of population dynamics. These 
equations are intended to describe the ecological mechanisms ruling the 
observed dynamics of a system, and radically excel formalized descrip- 
tions such as regression formulae, which feature widely in the modeling 
literature. including that for forest modeling. 

The various forms of forest (or stand) dynamics models can be 
subdivided. for our purposes, into three types (for a full review and 
details see Shugart (1984)): 

(a) Models with zero-level averaging over a tree population. They are 
based on the dynamics of individual trees with individual coordinates: 
spatial models in the exact sense of the word (Newnham & Smith 1964; 
Ek & Monserud 1974; Mitchell 1969; West 1987; Gurtzev & Korzukhin 
1988). 

(b) Models with an intermediate degree of averaging. They operate 
by considering individual trees but without their spatial location (Botkin, 
Janak & Wallis 1972; Shugart & West 1977; Shugart 1984; Huston & 
Smith 1987; Leemans & Prentice 1989), and form the family of so called 
'gap models' (see Shugart & Prentice, this volume, Chapter 12). 

(c) Models with the maximum degree of averaging. They operate by 
considering groups of trees-parts of a whole population — chosen depend- 
ing on the particular task (usual subdivisions are by size and/or by age) 
(Ek & Monserud 1979; Kapur 1982; Korzukhin, Sedych & Ter-Mikaelian 
1987, 1988; Tait 1988; Korzukhin, Antonovski & Matskiavichus 1989; 
Prentice er al. 1989; Kohyama 1989). 

Usually, different modeling techniques are more effective for solving 
one type of task and less effective for others. Approach (a) is suitable for 
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even-aged stands, for detailed investigation of stand structure, and for 
testing of different models of competition. As a rule, it demands much 
input information and is too bulky for parameter estimation by means of 
comparison of the model and empirical behavior. Intermediate case (b) is 
the most appropriate when the stand under consideration is multi-aged 
and composed of many (say, three or more) species. It is important to 
note that whereas models of types (a) and (b) are purely imitative, models 
of type (c) sometimes allow analytical approaches because they have a 
form of dynamical equations, 

Approach (c), a review and examples of which are given in this chapter. 
has been somewhat moved aside by the individual-tree technique. How- 
ever, we believe that this classical approach initiated by Lotka and 
Volterra, is complementary to those of (a) and (b), and has not exhausted 
its potential in theoretical and applicative tasks (for a contrasting view- 
point see Huston & Smith (1987) and Huston, DeAngelis & Post (1988). 

There are some fields of forest dynamics where dynamical equations 
are indispensable. For example, the dynamics of grass and moss popu- 
lations as part of the forest ecosystem cannot really be described at the 
individual-plant level; the same is true for all non-living components of 
forest ecosystems. such as dead organic material, nutrients, etc. In many 
other fields, dynamical equations need not be used but may be appropri- 
ate and helpful. For example, the inevitable stochastic element in 
population dynamics when an individual-tree approach is adopted leads 
to difficulties in parameter estimation: the majority of calculative 
methods demand determination of parameter derivatives that are almost 
impossible to find in a stochastic system or can be achieved only with low 
accuracy. Dynamical equations, however, are ideally suited for the 
calculation of the small trajectory modifications under small parameter 
variations that are needed for finding the derivatives. 

It is clear from these general modeling considerations that this tech- 
nique can be most successfully applied for comparatively simple popu- 
lation systems. The technique results in numerical solving, effective 
parameter estimation and — the most attractive theoretical aim — possible 
analytical results. Further, many boreal forest ecosystems are composed 
of a small number of tree and above-ground plant species, so the use of 
dynamical equations in this context can be of applied interest. In contrast, 
dynamical models for multi-species tropical and deciduous forests consist 
of many equations and parameters, so the deterministic behavior of tree 
number trajectories becomes so complex that it is near-stochastic. As a 
result, the advantages of a simple description by means of a dynamical 
equation technique over an individual-tree approach cease to exist. 
Besides the aim of the review and development of population level 
models in forest dynamics, we have attempted to illustrate how, in 
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our opinion, the approach to modeling should develop. At present, 
modeling appears to be more akin to art than science with each model 
being more or less specific, and each modeler choosing his own empirical 
rules of advancement. The view we advocate can be expressed as ‘Try to 
use simplicity before sinking into complexities’. For example, if the 
requirement is to analyze the behavior of the overall population-level 
characters only these, if possible, should be used as state variables, 
without size or age distributions, individuals, etc. If single-level, horizon- 
tal crowns give reasonable results, these should not be changed to 
vertically distributed ones. If a single-variable equation of trajectory of 
tree growth is required as output, use a single-variable equation of growth 
until it gives us what we want. The model should be made more complex 
only if it does not provide the required output. Many examples of 
unnecessarily complicated approaches are found in the forest and ecology 
modeling literature. 

In general, much modeling is performed on a ‘trial and error’ basis. 
Usually, the ecological mechanisms and the level of detail that is to be 
included in order to achieve the required information on a population or a 
forest system are not known. As yet, no attempts have been made to 
answer this basic question, although the whole issue is central to any 
modeling. Only when the necessary and sufficient relationships are 
known can one affirm that property A is explained by property B, that is, 
B is the reason for A's origin, Only a model constructed on this basis gives 
‘explanation’ and can produce accurate and reliable results; all others are 
incomplete to some degree. 

Forest modeling is currently in an active phase of its development. 
Both the expectations and the resources involved are high, and we 
believe that a certain degree of methodological accuracy should be sought 
in order to optimize output and to make the results more reliable. 


13.2. Major divisions of forest modeling 


The problems of the formulation of adequate forest models are in 
principle similar for individual-tree and dynamical equation techniques. 
In order to construct any model of a stand one should describe, for multi- 
aged populations: (1) the free growth of one tree; (2) the mechanism of 
competition between trees; (3) the dependence of the viability of one tree 
upon available resources and environmental factors: and (4) conditions 
necessary for seed viability. Of these points (1) and (2) are commonly 
required for individual-tree and dynamical equation techniques; (3) and 
(4) can be embodied in either deterministic or stochastic forms. thus 
giving the two discussed approaches to forest modeling. 
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Given these major divisions, in various forms and details, the modeler 
should be able, in principle, to pick out a model with desirable properties. 
However, at present we are some way from this ideal situation; there are 
still relatively few basic models and these have not been efficiently tested 
to allow them to be used in such a way. 


13.2.1. Tree growth 


The most physiologically correct and apparently most widespread 
approach to deriving tree growth equations involves the use of the 
concept of carbon balance (Richards 1959; Pienaar & Turnbull 1973; 
Aikman & Watkinson 1980; Bichele, Moldau & Ross 1980; Jarvis, 
Edwards & Talbot 1981; Waring & Schlesinger 1985; Mákela 1986; West 
1987; and many others). The simplest approach uses only one variable, m 
(total tree biomass), so that the balance takes the form 


m = S,(m)- A(R) — Re(m) — De(m). (13.1) 


Here m is the biomass growth rate, S, is total leaf area, A is specific 
assimilation rate. R is resources of photosynthesis, Re is respiration, and 
De is decay losses. 

In spite of the obviousness of Eqn (13.1), it is useful to note that when 
we consider long-term (e.g. one year) carbon balance, the difference 
S;A — Re cannot be less than zero. In fact, the balancing of Eqn (13.1) 
does not really represent the balance of tree biomass, but the balance of 
tree free assimilates (plus decay): free assimilates = 5;A — Re when 
S;A > Re and = 0 when S;A < Re. 

Since Eqn (13.1) uses only one variable for the growth description of a 
whole tree, 5, , Re and De are expressed as functions of m. It is common 
for more than one variable to be introduced (Agren & Axellson 1980; 
Kull & Kull 1984; Oja 1985; West 1987; Makela 1988; and others). In such 
cases, the authors are immediately confronted with the problem of 
assimilate distribution. Empirical difficulties are obvious, but theoretical 
approaches, which should be based on the optimization theory of growth, 
have not been satisfactorily developed. 

In this area, the work by Rachko (1979) is of particular note. In this 
work, the roots, stem and leaves of a tree are treated as independent 
variables, and the tree is assumed to be growing with the requirement of 
maximization of biomass increment during one year. The resulting 
description of photosynthesis is fairly realistic. Among recent works, we 
can mention Mákela (1988) who considered optimal tree height growth, 
and Korzukhin (1985) who introduced competition into the optimal 
growth task. The work by Korzukhin & Ter-Mikaelian (1987), which 


Thus, we believe that as the first step in a population level approach a 
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growth model description should take the form of a one-variable growth 
curve (e.g. D(t) — diameter or m(r)), and we recommend the use of the 
approach shown in Eqn (13.1) for this purpose. 

Models differ in the forms of A= A(R), Re- Re(m) and 
De — De(m). Here we shall not consider the vast literature on corre- 
sponding models, but note only that the simplest and most widespread is 
the hyperbolic function A(R;), where R; is intensity of light. However, 
one needs to take care in the application of this classical light curve to a 
whole tree growing over a long period: the mechanism of photosynthetic 
adaptation can make the function A(R;) almost linear (Tcelniker 1978). 
Respiration of growth is proportional to 5,.4 and can be taken into 
consideration by a coefficient before the term 5$; A in Eqn (13.1). In the 
simplest case. other components of respiration are taken to be pro- 
portional to m (Bertalanffy 1957; Richards 1959; Pienaar & Turnbull 
1973; Aikman & Watkinson 1980; and others). In more complex cases, 
some elements of the respiration structure are introduced (Bichele, 
. Moldau & Ross 1980; Mákela 1986). De(m) is usually taken as being 

proportional to m. s 

Some authors (Botkin. Janak & Wallis 1972; Reed 1980; Shugart 1984; 
Leemans & Prentice 1989; and others) have used a semi-empirical 
analogue of the balance equation (13.1): 


m= SA - (1 — mim™). (13.2) 


Owing to Eqn (13.1) the value of m™** (maximum tree biomass (the root 
of m = S, (m)A(R) — Re(m) — De(m) = 0)) depends on the assimilation 
rate A. For example, A and m"** are less in poor conditions. On the other 
hand, owing to Eqn (13.2), m™* is independent of A, and only the 
growth rate depends on A. 


13.2.2. Tree viability 


The problem of tree viability modeling consists (a) of defining the 
concrete form of function 


Y(x), (13.3) 


that is, the probability of survival of one tree over, one may suppose, one 
year, and (b) of defining the list of arguments, x, of the function. 
Following the usual subdivision (Semevsky & Semenov 1982; Tait 
1988), viability factors are divided into two large groups. The first of 
them, denoted as w, is connected with stochastical, often catastrophic, 
natural disturbances like fire and drought (Waring & Schlesinger 1985) so 
w — (frequencies and magnitudes of disturbances, climatic parameters, 
etc.). These factors are the arguments of the density-independent com- 
ponent of viability, W(w). The second group. denoted as v, arises from 


Population-level models of forest dynamics 339 


Viv) 


1 


D, Susa: ADS WG 0 v lv 


(a) (b) 

Fig. 13.1. Qualitative behavior and arguments of viability function Y 
for a single tree. (a) Empirical (1) and model (1), (2) forms following 
field data and models (see text for references and denotations). W is 
viability of mature, maximum-defensed tree. (b) Qualitative behav- 
ior of density-dependent viability component, V. on argument v — 
AIA maxi V is recommended for Base model 1 (83.3) and used in Base 
model 2 (§ 4.3): vy, = Amin/Amax: where Amin is maintenance 
respiration of leaves or a whole tree. 


tree-tree permanent interrelations in the community. and gives argu- 
ments of the density-dependent component of viability, V(v), where 
v — (all tree state variables and density). So we can represent x as 
X = (W, V) and write 

Y(w, v) = W(w) - V(v). (13.4) 
Below, we shall discuss the problem of V(v) dependence only, which has 
been investigated much less than W(w). 

Several attempts have been made to establish the form of V on the 
basis of experimental data. In most cases, arguments of V used have been 
dimensional characteristics of the tree, its growth rate, and a parameter 
describing the pressure of competition. The dependence of V on these 
arguments is usually quite marked (Moser 1972; Monserud 1976; Glover 
& Hool 1979; Hamilton 1980; Buchman, Pederson & Walters 1983; 
Blagovidov 1984) as may be seen in Fig. 13.1(a). Glover & Hool (1979), 
Hamilton (1980), and Buchman, Pederson & Walters (1983) used the 
exponential function 


V(v) = V[1 + e^ *(?], 
where v is a set of arguments, different for each of the cited papers 


(v = D, G, D, D/D, competition index), and g a linear function; here 
and below D is average tree diameter, Nis tree number, and G = (D?N)/4 
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is the stand cross-sectional area. Moser (1972) found from empirical 
observations that V = 1 — exp (—aD), and Kohayma (1989) found that 

— ] — aG. For a description of empirical data, K. J. Mitchell (1969) 
used the argument v = S./SP**, where S, is the real crown projection 
and S™* the same for a free-growing tree. This argument means that a 
suppressed tree ‘remembers’ its possible size if its growth is no longer 
suppressed. The same effect is proposed by Ek & Monserud (1974). 

So far. the most profound approach to the problem has been the 
theoretical construction of Korzukhin & Ter-Mikaelian (1987) (see also 
Antonovski, Korzukhin & Ter-Mikaelian (1987)). Tree life function 
‘viability’ is considered as being one of a number of activities (growth 
and proliferation in these models). The tree is assumed to be a system 
which distributes resources (free assimilates in these models) among 
the activities, and of high priority is the maximization of total seed 
production. This is proportional to accumulated leaf area and tree 
longevity, which in turn are determined by the amount of assimilates 
directed to defense. This distribution task was formalized and its sol- 
ution showed that the share of assimilates directed to defense decreases 
monotonically with age. A trajectory of optimal growth m(t) and realis- 
tic age behavior of seed production was also obtained: this was equal to 
zero until a certain age and thereafter increased monotonically. The 
results of such an analysis for a free-growing tree and for a population 
of even-aged trees were qualitatively close. Despite being sound from 
the biological point of view, this approach does not lend itself to 
practical applications. 

More practical modeling realizations of empirical observations can be 
divided into discrete and smooth, the first type being appropriate to 
individual-tree level models and the second to all kinds of models. The 
preferred argument of V is usually D because of its relatively easy field 
measurement. 

The application of a step-like function V was initiated by Newnham & 
Smith (1964) who used V = 1 when D > const and V = 0 otherwise. This 
was used by Aikman & Watkinson (1980) taking $ instead of D. Mákela 
& Hari (1986) took V = 1 when $> 0 and V = 1 — const $,/5, other- 
wise. In a gap-modeling technique (Botkin, Janak & Wallis 1972; Shugart 
1984; and many others), V = 10 77 when D > 1 mm yr^! and V = 10 9? 
otherwise, where T is the lifespan of the species. K. J. Mitchell (1975) has 
taken V = 1 when NH? < const and V = value sustaining NH? = const 
otherwise (H is average tree height). Tait's (1988) model is the only one 
that uses density and growth rate simultaneously as viability arguments: 
V = const — NS,. In West's (1987) model a tree dies when total assimi- 
lation becomes less than respiration of the leaves, S; A — rj m, < 0, that 
is equivalent to A < const, because S; = mz. 
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The smoothed forms of V are also varied. Korzukhin, Sedych & Ter- 
Mikaelian (1987, 1988) applied in their quantitative model . 


V =(R,/RP™)?, 


where R; is the light available for a tree, and a very similar argument was 
used by K. J. Mitchell (1969) (cited above). Zemlis & Shvirta (1987) took 
V—1- GN, and Cherkashin (1980) took V = v* exp (—1/v), where 
v = D/D. In their population-level model (see $13.4.3), Antonovski et al. 
(1989) and Korzukhin, Antonovski & Matskiavichus (1989) used 
V= 1 — const (1 — v), where v = A/Amax is the relative assimilation 
rate of a tree (Fig. 13.15). Pukkala (1987) applied V = 1 — G/tasseedling 
viability. The first case of a smooth viability curve V(D) in gap modeling 
was presented by Leemans & Prentice (1989). 

To conclude this short review, we note that the way forward as far as 
definition of viability is concerned must lie in the use of evolutionary 
optimization principles (Pianka 1978; Semevsky & Semenov 1982). As 
yet, their application has not led to algorithms for the definition of V 
which are reliable enough to be used in a forest dynamics model (see 
Korzukhin & Ter-Mikaelian (1987) and Antonovski, Korzukhin & Ter- 
Mikaelian (1987)). 


13.2.3. Competition 


Combined resource utilization leads to interaction between trees and 
unifies individuals in a stand, so competition can be considered as one of 
the two central phenomena that determine the dynamical behavior of a 
forest (the other is interaction connected with reproduction). 

In essence, the problem consists of defining tree density-dependence of 
two functions: assimilation rate in Eqn (13.1), A — A[R(N)], and viabi- 
lity in Eqn (13.3), V = V[v(N)]. Density can be the average number of 
trees per unit area, in population-level and local approaches, or several 
neighboring trees around the central tree in the individual-level 
approach. 

As in the cases of tree growth and viability, there are two ways of 
introducing density, namely formal and mechanistic approaches. The 
first is widely represented in forestry modeling and involves the formal 
introduction of ‘indexes’ of competition, CJ. When the individual-tree 
level is considered, C7 are artificial functions of diameter, height, cross- 
sectional area, and biomass increments of neighbors of the tree in 
question (Newnham & Smith 1964; K. J. Mitchell 1969; Ek & Monserud 
1974; Diggle 1976; Britton 1982; Smith & Bell 1983; Cennel, Rothery & 
Ford 1984) or of local density, that is number of trees (Plotnikov 1979), 
or of quite formal parameters of the central tree or its neighbors (Laessle 
1965; Fires 1974; Lorimer 1983). When CI have been chosen, the 
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correlation between the diameter increment of the central tree and CJ is 
usually calculated, with the aim of testing the competition model. 

When description is at the population level, that is dynamical 
equations are used, CJ are some function of mean population variables 
such as D/D (Cherkashin 1980), GN (Zemlis & Shvirta 1987), H^N (K.J. 
Mitchell 1975), N (Chjan & Chjao 1985), e, (Tait 1988), or G (Kohyama 
1989). 

The above approach is rather formal: in particular, C/ cannot be used 
as an argument of A in Eqn (13.1) if one wants to interpret A in 
physiological terms. 

The other, mechanistic, approach deals with the physical resources. R, 
available to a tree and tries to calculate their dependence upon density 
(taken again in two senses), that is to define the form of function 


R(N, r), (13.5) 


where r are some dimensional characteristics of neighboring trees. or 
trees of the whole population. From general considerations it is clear that 
the results will depend on the type of resource. tree morphology. and 
spatial pattern of the trees. 

Our understanding of competition for light is considerably ereater 
than that for soil resources. Usually, a uniform distribution of tree foliage 
is assumed (a representation of the producing layer) that immediately 
leads to the application of the Lambert-Beer penetration law (Botkin. 
Janak & Wallis 1972, plus all gap modelers; also Ross 1975: Mákela & 
Hari 1986; Oker-Blóm 1986). This law states that the amount of light at 
height level h is given by 


RA) = Ryo exp [—yo*(h)]. (13.6) 


where y is the extinction coefficient after transmission through one leaf 
layer, which can depend on h, o^ (A) denotes the amount of foliage layer 
above level h, and Rp is the initial light flux. 

However, the sancta simplicita of Eqn (13.6) is lost after the first 
refinement of the uniform model. Foliage in a stand is organized in tree 
crowns. Let us consider the question of light available to an average tree 
in a population (Korzukhin & Ter-Mikaelian 1982). Each tree has a 
horizontal monolayer crown with area S; , and n(A) is the number of trees 
with height A. It is important that stem bases are distributed randomly 
over the plane, a requirement of Poisson's law. The extinction coefficient 
after one foliage screen transmission is equal to y. By means of a 
geometrical probability technique, it was shown that under these con- 
ditions the classical formula (13.6) is true, where o^ (h) = f; S,n(z) dz. 
The lack of a Poisson distribution, however, violates Eqn (13.6): this law 
is necessary for the fulfillment of Eqn (13.6). 
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Another example is provided by the same system but with, say, two- 
layered crowns. Each layer is close enough to the other to avoid 
intersections with other crowns and has area S., so S} = 2S,. Thus, the 
whole-tree absorption coefficient becomes equal to [= 1 — (1— y). 
Following our approach, R,(/) is given again by (13.6) but in the form 

Rh) = Rp exp [Taz (A)], 
where o; is the area of crown projections above level A, which in our case 
consists of 0° (A). Finally, we obtain 
Rh) = Ryo exp [7316 (h)], 

that is the analog of Eqn (13.6) with effective extinction coefficient 
y -[1 — (1 — yy V2. a new violation of the uniform distribution model. 
Obviously k-layered crowns will give y'= (1— y)*]/k. The physical 
reason for this violation is clear: being organized in crowns, foliage 
screens have lost their freedom to move independently of each other, and 
only the whole crown (k screens together) has reserved this ability. 

Several such grouping foliage models have been developed for crowns 
and stands (for review and references see Oker-Blóm (1986)). If the 
possibility of crown intersections is to be taken into account, the analyti- 
cal calculation of R,(/1) becomes almost impossible. 

Although root competition is just as important in forest dynamics as 
competition for light, no one model for the former has been developed 
that corresponds to that for light. Most approaches do not use soil 
resources as explicit variables, although some do. As a rule, authors 
determine some ‘qualitatively true’ functions of total biomasses of 
competing trees (moss and grass if needed), represented as B(M); it is 
usually assumed that root competition has a suppressive effect, so dB/ 
dM « 0. These functions are then used as multipliers to assimilation rate 
function A (e.g. McMurtrie & Wolf (1983) and all gap modelers). The 
most recent example of this (Shugart 1984) is a multiplier for the 
individual-tree growth equation 

B(M) = 1— M/M™*, 
where M™ is determined from field observations. 

We can offer here a more realistic approach to the modeling of root 
competition, based once again on geometrical probabilities (Korzukhin 
1986). Let us consider a population of identical screens, each having a 
thin, plane root system with area Sg. All are disposed at one level as is 
often the case for the boreal forest zone, so the process of competition is 
two-dimensional. Root systems are distributed according to Poisson's law 
over the plane, and from a unit area of the medium an amount Ro of 
resources (water, oxygen, etc.) is available. These resources are divided 
equally among all root systems which overlap at a given point; that is, 
from area Sp. which is not overlapped by other roots, the tree obtains an 
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amount Rs; = Sai Rso of resources; from area Sp>, which is overlapped by 
one neighbor. it obtains Rey = SgoRso/2, etc. So, the total resource 
available to a particular tree is 


p = Sg, > Sri = Sg. 


i=l 


We are interested in the average amount of resource 


Rs = E W;R so. 
i-i 
w; are average areas of the ith overlapping points. It can be shown that 
under Poisson's law Sg/Sg =å" e^ ^/(i — 1)!, where À = SRN. that is 
average coverage of unit area (root area index). Finally we have 


&,-T2p — exp (—SaN)]. : (13.7) 


Magnitude R can now be used as an argument in the assimilation rate 
function. together with light: A = A(R;. Rs). One variant of the function 
was developed by Gurtzev & Korzukhin (1988) and applied in an 
individual-tree model of a linear pine stand: it was shown that taking into 
account the process of root competition improved appreciably the quality 
of growth description. 


13.3. Even-aged monospecific stands 


13.3.1. Introduction 


If initially we do not take into consideration the problems connected with 
the migration and establishment of seeds, even-aged stand modeling 
consists of the same basic elements as modeling of multi-aged and multi- 
specific stands. In both cases one must be able to formalize the processes 
of competition of individual tree growth (including, perhaps, changes in 
morphology) and of tree mortality. However, the field data relating to 
even-aged stands are much more accurate and numerous than those for 
multi-aged ones. 

So, it seems reasonable to adjust a population-level model on the basis 
of even-aged behavior first of all, and if this is successful, progress to 
multi-aged behavior. Surprisingly, this apparently obvious way of model 
development has not been accomplished until now; even- and multi-aged 
stand models have developed independently. 


13.3.2. Empirical behavior of even-aged stands 


The required data can be taken from numerous observations at perma- 
nent plots. We shall restrict ourselves to the level of description of a stand 
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which uses only the average characteristics of a tree, but not tree 
distributions. In this case, the system is described by the following 
variables: tree biomass m, diameter D, height H, leaf area 5; , and seed 
production p (per year). The population variables are: tree number N, 
total biomass M = mN, leaf area o = S, N (or close to o is cross-sectional 
area G), and total seed production P — pN. 

For the purposes of model development it is necessary to know the 
behavior of these variables under variation of initial density N(0) and 
ecological (site quality) parameters g. A brief summary of this behavior is 
given in Korzukhin (1986) and in McFadden & Oliver (1988). 


Time behavior 


Individual-tree variables have a simple form of monotonic sigmoidal 
functions: population variables (except N) are nonmonotonic and go 
through a maximum (Fig. 13.2). After crown closure, a specific system 
invariant arises, connected with the maximum amount of leaves, 9^?*, 
which can be achieved in given ecological conditions: 


olt) = S, (r)N(rO) = o™. (13.8) 


Density behavior 


Individual-tree variables m, p and D[t. N(0)] taken at any given moment 
decrease monotonically with the increase in initial density N(0) (although 
H{t, N(0)] is sometimes not monotonic). Maximum differences are ob- 
served somewhere in the middle of the set of trajectories for m (Fig. 
13.3). 

If one considers tree number dynamics under N(0), one can observe 
the effect of ‘forgetting’ of initial conditions (Fig. 13.3), that is 


N(0) 


8 t 


Fig. 13.2. Qualitative dynamics of major variables of even-aged 
systems; curves with the same type of behavior are unified; p = m N?? 
(see text for definitions). 
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p! 


Fig. 13.3. A set of even-aged svstem trajectories obtained under 
various initial densities (see text for definitions). 


N(t. NS(0)/N(t. N,(0)] — 1 for all N, .(0). This “drawing together’ effect 
means that tree mortality depends upon density (and increases as density 
goes up). 

Population variables M. P. G and o(r. N(0)] increase at first with 
increased N(0). Owing to somewhat sparse data (Results of Experimental 

. 1964; Hirano & Kira 1965; Buzykin 1970; Redko 1978). their 
maximum values M™**, etc. begin to decrease under very large N(0) (Fig. 
13.4). We call this effect ‘overcrowding’. The time at which these 
variables reach their maximum values, 47°G ,[N(0)]. decreases with 
increasing N(0) (Fig. 13.4). 


Site quality 

These effects are obvious for individual trees. At the population level, 
tree numbers are smaller in the better site conditions (under the same 
initial conditions) (Fig. 13.5). In spite of this effect, total biomass and 
other population variables increase. 


Tolerance effects 


Under the same initial conditions, tree number is greater for more 
shade-tolerant species (Fig. 13.5). 


13.3.3. Dynamical models of even-aged stands 


Surprisingly, there are no even-aged stand models which are able to 
describe all properties of the system enumerated in Section 13.2.2; there 
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(a) 


Mit ha^!) 


(b) 


G(m? ha) 


t(vears) 


Fig. 13.4. Trajectories of population biomass M (a), and cross 
sectional area G (b), (which is proportional to c) under initial 
densities N,(0) = 1000 (pluses), N3(0) = 2400 (filled circles) and 
N3(0) = 10000 (open circles); data from permanent plots with Pinus 
sylvestris in Moscow region (Results of Experimental ... 1964). 
"Overcrowding' effect consists of MP?* <MP® under N;(0) > N;(0), 
whereas all G™ are approximately equal. so only r@** changes. 
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N,(0) 


p> 
wN 


t 


Fig. 13.5. Given the equal initial numbers of trees, N, (0), N; (t) falls 
rapidly for ameliorated site conditions (1 — 2 — 3) and with decreas- 
ing shade tolerance (1' — 2' — 3') (qualitative behavior). 


exist only some formalized constructions relating to different parts of the 
whole picture. 

A number of works are traditionally devoted to tree number behavior 
only. Among them, there are linear equations N = —d(t)N where 
d(t) = 1 — V(r) is specific mortality as a given function of age (Terskov & 
Terskova 1980; and others); obviously, this equation does not ‘forget’ 
initial conditions N(0) and consequently cannot describe the ‘drawing 
together’ of trajectories. In addition, there are nonlinear “Volterra-type’ 
equations, e.g. N= —aN[1 — (b/N)*] (Chjan & Chjao 1985). 

On the other side, there are ‘productivity’ models, operating with one 
variable M: M = aM“ — bm (e.g. Pienaar & Turnbull 1973; Budyko 
1977; and others); obviously, this equation cannot describe M(t) going 
through its maximum. 

The majority of modeling efforts in this area have been devoted to a 
special class of model based on different forms of the so-called ‘3/2 law’ 
invariants (Reineke 1933; White & Harper 1970; Hozumi 1980; Lonsdale 
& Watkinson 1983; Zeide 1987; and many others): 

p = yN" = const. (13.9) 
where y = S}, m, D, or H. It is likely that the possible population level 
‘invariant’ is interrelation (13.8), so Eqn (13.9) is, in fact, Eqn (13.8) 
rewritten in terms of other variables. 

If our aim is to model the whole scenario adduced in Section 13.3.2. 
interrelation (13.9) appears to be very restricted. Firstly, it cannot be 
related to the part of the trajectory M(t) which lies behind M™™. If we put 
y=m (classical variant of Yoda er al. (1963)), then p= mN? = 
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MN?-! = const; because Ñ < 0 this leads to M > 0, because (a — 1) is 
greater than 0, so that Eqn (13.9) can be true only before M™**, where 
M = 0. Secondly, it cannot take into consideration the important depen- 
dence of the trajectories on initial density N(0), as this value is not 
‘remembered’. 

The most developed approach of Hozumi (1980) uses, besides in- 
variant (13.9), non-autonomous dynamical equations for individual-tree 
growth, and looks artificial. 

All types of dynamical analysis of even-aged stands based on the 3/2 
law give us a ‘semi-model’, which is intended to substitute the formu- 
lation and analysis of a complete nonlinear system of dynamical 
equations aimed at modeling the whole trajectory of the system. It seems 
to us that in this area there is a certain ‘magic’ attached to simple 
formulae, and that researchers have now extracted from them all that was 
possible. 


13.3.4. Base model 1 of even-aged stands 


Below we use the simplest model elements of growth, viabilitv and 
competition for the composition of even-aged model 1. The aim of this 
model is to describe the maximum dynamical properties of even-aged 
stands enumerated above. 

Calculations using growth equation (13.1) for free-growing trees and 
trees experiencing competiton show that it works well when one simply 
takes 


Si (m) = am?, (13.10) 
and when 
Re(m) + De(m) = cm (13.11) 


(see, for example, papers in Richards 1959; Fires 1974; Gurtzev & 
Korzukhin 1988). Interrelation (13.10) can be derived from two well- 
established ties: S; — DÊ, 1.5 < 8 - 2 (e.g. Mohler, Marks & Sprugel 
1978; Mákela 1986), and m — D?, 2.5 € 3 (numerous yield tables). As a 
result we obtain Eqn (13.10) with 0.5 =a x 0.8. 

Following (13.4), total viability Y — WV. In the construction of V, we 
are working from the basis of the following propositions: (1) it should be 
based on a resources approach; arguments such as D are considered to be 
indirect reflexes of the real viability mechanism, which is in essence the 
use of tree resources for defense and repair; (2) we will abandon the 
‘memory’ arguments such as m(t)/m""*(r) for the sake of simplicity, since 
they demand a second growth variable (for m™**(t) for example); (3) 
because photosynthesis is the central process of tree resource production 
we suppose that assimilation rate A should be the argument of V in Eqn 
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(13.4) (normalized for suitability): 
v = A/A max 


Following (13.1) and (13.5), A depends on N through R, that is 
A = A[R(N, m); when N = 0, R = R™™, A = Ajax. V = 1, and V = 1, 
we obtain a free-growing tree. 

Finally, we have base model 1 (specific death rate is equal to 1 — Y): 


i = {1 — WV[A(R(N, m))/Amax]}N. 


j (13.12) 
m= S,(m)A[R(N, m)] — em. 


where S, (m) is given by Eqn (13.10), W < 1, V < 1, A is given by any 
model of photosynthesis, and R is given by the model of competition. 
Below we adduce a simplified analytical example of the use of Eqn 
(13.12). Again, we consider the competition for vertically directed light 
among a population of Poisson-distributed trees with horizontal crowns, 
thinly spread over height. This gives (Korzukhin & Ter-Mikaelian 1982): 
Rj(N. m) = Ryo exp [75S, (m)N]. (13.13) 
Consider the situation when competitive interaction is weak. that is 
magnitude x = 7/5; (m)N «& 1, and undertake Taylor's expansion of A 
near x = 0 and V near v = 1: 


A(x) = A(0) — (94 


X—..4 
dx 


Vv)» V1) - =) - ..., 


where A(0) = Ajax, V(1) = 1 and the argument of V is 
1 |da 
Zins 


dx 


Assuming that all first derivatives are not equal to zero, |dA/dx| = a, = 0, 
|dV/dv| b; = 0, we obtain a system of first approximation 


v= A(X)/Amax = 1 — 


N= |1- w(i- ab, 1 sw. 
max 


(13.14) 
ay 


m = S,(m)A max (1 = x|- em, 


*^max 

Now suppose that density-independent mortality is absent, W = 1. 
respiration and decay losses in m are negligible, and a4x/A max < 1. that is 
growth is free. The result is an idealized system 


N2-m'"N; m= qm^ 
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(r = aajbilAq4.; ¢ = @Amax) and its solution under initial conditions 
N(0), m(0) = mo: 


m(t) = [(1 — a)qt + my ^]"^*; (13.152) 
1 X = —a ie 
N(t) = Is: - = + A [(1 — a)qt + mj-*]! | .. (13.15b) 


It is clear that solution (13.15a) quickly forgets initial condition mo (mass 
of seed), so the formulae are simplified: 


-1 
m(t) = [1 — ag] *, NO) = I 4e [rt ay (13.16) 
No q 
We also need a generalized population variable of the form F = m^N. 
which, from (13.16), is equal to 
pj" a 


Fy = 13.17 
OTN, +r (RN 


(ne ((1 T a)q]^^-*, rs = = r(1 = qe agra), 


F represents M, g and G when f is taken as needed. and under B < 1 it 
goes through a maximum when ¢ = r£^*: 


Faux = No PAL AE (13.18) 
L- uve " 

mp = —— Ana fa(é) (13.19) 
Ni 


where we have included in f, and f» all dependencies upon the rest of the 
parameters, which are not of interest here: 5 = (a. a, b4, y). 

Equations (13.15)-(13.19) correctly describe many properties of an 
even-aged system: 

]. From (13.155), density N[t, N(0)] forgets initial conditions N(0) 
and converges to a ‘magistral’ trajectory (Fig. 13.3): 


1 1 
abya 2 A2a- Vi-a dia 


2. When site conditions ameliorate, A,,,, increases and N(/) dimi- 
nishes at any given í (Fig. 13.5) if one supposes a > 0.5 (a fairly realistic 
condition). 

3. If we consider more shade-intolerant species, then the reaction of V 
to diminishing of A becomes sharper, that is b, (and consequently 7) 
increases and leads again to N(t) diminishing (Fig. 13.5). 

4. When No increases, Fmax increases also (for f < 1), that corre- 
sponds to the behavior of o™*(No and G™*(No) (Fig. 13.4) under a 
relatively small N(0). 


N(t) = 
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5. Magnitude r™*(No) (Fig. 13.4) also decreases under better site 
conditions (when a > 0.5). 

Taking into account respiration and decay in the tree growth equation, 
m — qm^ — cm will make M(t) go through a maximum because m(t) 
becomes finite but N(r) — 0, and addition of competition will obviously 
give a bunch of trajectories m[r, N(0)]. 

Finally. the only dynamical effect that cannot be given by base model 1 
(13.12) is the ‘overcrowding’ (Fig. 13.4) under large N(0); it can be shown 
that this effect demands at least three dynamical variables. 

One can easily see that the form of system (13.14) and all results do 
not, in fact, necessarily depend on the resource under competition: they 
are true both for light, and for soil resources. Bearing in mind the 
proposition about weak interaction, one can take Taylor's expansion of 
Eqn (13.7) and repeat the same calculations with x = gSp(m) << 1. 

The topic of numerical applications of system (13.12) is worthy of a 
separate paper. It seems to us that the above considerations show clearly 
the ability of base model 1 to be used for natural forest modeling. 


13.4. Multi-aged stands 


In the array of successively complicated elements in the field of forest 
dynamics, multi-aged monospecific stands are the simplest elements, 
ones which can be related to real natural forest. The only process that 
needs to be added here in comparison with an even-aged population is the 
origin of seedlings. This gives us a usual population demographical 
system with a complete collection of dynamical processes. 

Itis necessary, then, to undertake a model description of the combined 
dynamics of the set of age cohorts which represent the whole population. 
There are some variants of the mathematical embodiment of the popu- 
lation which are distinguished by age and consequently by tree size and 
we shall review them briefly. 


13.4.1. Construction of demographical models 


There are a number of similar ways to formalize a population's age 
dynamics. The first approach was illustrated by Von Foerster (1959) who 
considered n(¢, t) (the number of individuals with given age r at time t) 
and processes of birth and death, ignoring growth. This corresponds to 
the model 


ms c on = —d(t.r)n. n(t, 0) = f p(t. z)n(r. z) dz. (13.20) 
ðt aT 
Here d and p are specific death and birth rates; note that they are 
independent of population density n, that is the model is linear. This 
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simple approach has now been exhausted from the mathematical point of 
view (Sinko & Streifer 1967; Poluectov 1974) and is of no serious interest 
for forest dynamics. Extensions of model (13.20) have been developed in 
various directions. Gurtin & MacCamy (1979) introduced density- 
dependent mortality d = d(t, N) where N = f n dr and obtained analyti- 
cal results for partial cases of d(N). Sinko & Streifer (1967) considered a 
two-dimensional system, which operates by combined age and size 
distribution n = n(t, t, m) where m is any quantitative characteristic of 
the individual. Behavior of n in the most general case is derived by 
Kolmogorov’s or Fokker-Plank's equation 

on | dn, ð AA 5 

* ES = - — [g(t. T. m. n) n] d(t, t, m, n)n (13.21) 
with corresponding boundary conditions; here g — mis the growth rate of 
the organism. Competitive and other density-dependent aspects are 
taken into account by means of argument n in g and d. In this case, Eqn 
(13.21) is very complex for analytical consideration; Sinko & Streifer 
(1967) have examined a partial case when Eqn (13.21) is linear with 
respect to n and m, that is g = g(t, t) and d = d(t. 7). 

Another and more popular approach was proposoed by Leslie (1945), 
whose well-known matrix technique is the discrete analogue of the 
continuous-time model in Eqn (13.20) and is more suitable for solving by 
computer. Leslie dealt with age distribution n(t. 7). It is convenient for us 
to write out his model in a ‘cohort’ form and with a generalized variant: 


ae lrt1)5Y(tr.(Dnmtr), t= lese T= 1 (13.22) 


n(t, 1) = p(t, v, (n))n(t. x) 


where t and t are discrete, Y and p are viability and birth rate, and (n) is a 
generalized vector argument, (n) — n(t, 1), n(t, 2)... .. n(t, T) which 
describes density effects. 

Leslie (1945) and many others have used this model with 
Y.p = Y, p(t) only, that is, the simplest linear variant which enables 
powerful matrix analysis. For example, Enright & Ogden (1979) pre- 
sented detailed numerical analysis of one- and multi-species rain forest 
population systems by means of a linear case of model (13.22) with 
empirically determined birth and death rates. It is clear that in the case of 
forest dynamics linear equations can be true either for low densities or for 
short time periods until argument (n) = const. 

Among recent works we can note a two-species age-distributed model, 
an application of system (13.22), presented in Korzukhin, Sedych & Ter- 
Mikaelian (1987, 1988). This was applied to 200-year post-fire succes- 
sional dynamics in West Siberia. Growth curves for both species (birch 
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and Siberian pine) were fixed, so only age number dynamics were 
analyzed, that is, behavior of magnitudes nj(r, t) and nə(t, 7). Crowns 
were horizontal and light competition interaction was directed from 
higher to lower trees only, as in the gap-model approach. The dvnamics 
of two age packages (‘waves’) observed in the field were observed 
numerically (see also $13.4.2.). 

In forest dynamics applications the Lefkovitch (1965) approach is 
more popular; this offered the same technique for size distribution 
analysis (sizes of trees are measured much more easily than their ages). If 
one breaks the size axis m into Q intervals 


(Mo, 1). Ma). - - Q1 kto) (13.23) 
and takes all trees whose sizes belong to interval i. that is uj, € m € uj. 
then the system dynamics will be represented by the following scheme: 


zm) M nt, 2) x (13.24) 
Bich Pe» jenes ete 
I n(t + 1, 1) nt + 1.2) 


where n(r. i) is tree number in the ith size interval. g, is growth and d; is 
death rate. It is proposed, for the sake of simplicity. that changes in sizes 
are small and occur only between neighboring classes. Dynamical 
equations, being discrete analogs of differential equations (13.21). are 
easily written from the balance scheme (13.24): 


a(t + 1.i) = niti) + g(t.i— 1, (n))n(t. i — 1) 
= g(t, i. (n))n(r. i) — d(r i. $))n(r. i). i22.....Q: 
(13.23) 
n(t + 1,1) 2 n(r. 1) + Y pli, i (n)n(t. i) — gi. 1. ())n(r. 1) 


— d(t, 1. (n))n(t, 1). 


Here p is size-specific fecundity. argument (n) is analogous to that used in 
Eqn (13.22), and dependence of g, d, and p upon i means their depen- 
dence on size. In these equations, the magnitudes of size intervals Axi are 
considered to be included in functions g. d. and p. 

The most frequently used variant of the highly generalized model 
(13.25) consists of taking birth. growth and death rates as depending only 
upon size (class number i). Buongiorna & Bruce (1980) have applied this 
model to the task of productivity maximization in this linear variant. 
when g = g(i). Hartshorn (1975) has used this type of model for two 
tropical tree species, and Dyrenkov & Gorovava (1980) for one spruce 
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species. The value and restrictions of the linear approach were noted 
above. Kapur (1982) undertook an important extension of the task by 
introducing g = g(N), that is, a rough reflection of density-dependent 
factors. The aim was again to maximize total stand productivity. A 
nonlinear model by Ek & Monserud (1979) treated the same object — a 
six-species forest in Wisconsin — as was analyzed by their individual-tree 
model (Ek & Monserud 1974). This is the only modeling work where such 
an interesting comparison was made. Unfortunately, the authors have 
not adduced the explicit form of their dynamical equations. A complex 
nonlinear model was formulated and applied to a seven-species boreal 
forest in South Siberia by Cherkashin (1980). Recent work of Kohyama 
(1989) is devoted to the analysis of tree diameter distributions dynamics, 
n(D, t), in a rain forest population system which was close to an even- 
aged one; the model is presented by nonlinear system (13.25) with 
empirically determined death and growth rates as functions of G. 

A major weakness in the described matrix and continuous-time multi- 
aged models is the rather formal realization of individual-tree growth and 
competition mechanisms. From this point of view. individual-tree models 
are better developed. 

A third way to investigate age and size dvnamics is by means of a 
combination of Leslie's (1945) and Lefkovitch's (1965) matrix technique. 
which can be made by means of straight generalization of our even-aged 
base model 1 in Eqn (13.12). In other words, a multi-aged model can be 
constructed simply as the sum of one-cohort equations (13.12). This 
method was proposed and realized by Antonovski, Korzukhin & Mats- 
kiavichus (1989) and Korzukhin. Antonovski & Matskiavichus (1989). 
This gives the system to which we shall refer as base model 2: 


n(t + 1, 1) = Ep(t, t, (n), (m))n(t, T) + f, 
n(t + Lt + 1)= Y(T; v, (n), (m))n(t, x), (13.26) 


m(t+1,t + 1)= m(t. x) + g(t, v, (n), (m)), | ol rn a E 


Boundary conditions for the growth equation are m(t, 1) = m,. Variable 
m, indeed, can be substituted by arbitrary size characteristics of a tree. 
Argument (m) has the same meaning as (n). Note that Leslie's system 
(Eqn (13.22)) is the partial case of Eqn (13.26); in order to obtain Eqn 
(13.22), one would simply not take into consideration all equations for 
m(t, x), that is, one would suppose that growth curve m(t) is set. 

We have added one new element in base model 2 compared with Eqns 
(13.22) and (13.25), namely addendum fin the equation for first age class. 
This element describes seed influx in the system due to the background of 
seed that usually exists in forested areas and introduces a spatial aspect to 
the system analysis which is absent in the standard form of Leslie's and 
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Lefkovitch's models. It is interesting to note that gap models, on the 
contrary, use only the background of seeds and neglect seed production 
inside the considered area itself. 

Let us briefly discuss the comparative potential of size and improved 
age cohort approaches, that is, systems (13.25) and (13.26). Major 
differences are the following. System (13.26) is purely deterministic, 
giving one trajectory of biomass and tree number for a given cohort, so 
that this even-aged sub-system at each moment in time is described by 
only two magnitudes: m(t + k, k) and n(t + k, k); k is the age of the 
cohort born at moment r. In contrast, Eqn (13.25) describes spreading of 
sizes (and numbers, correspondingly) about the m-axis, even for one 
cohort. In reality, after k steps each cohort will be partly presented in all 
size classes from the first to (k + 1)th. So certain stochastic mechanisms 
are contained in Eqn (13.25); in its strictest sense, as it deals with size 
distribution, it is analogous to a Markov chain. 

The two systems also have similarities. Both give us size distribution, 
model (13.25) by the definition of n(t, i) and model (13.26) after simple 
summation over a given size interval (Eqn (13.23)): 


n(t,i)= Y n( v), TES TS Tq, 


T 


where rz; and v;,; are the ages whose tree size m; belongs to the ith 
interval, 4; € m; « u;.,. Both systems also give us age distribution. 
model (13.26) by the definition of n(r, t) and model (13.25) by means of a 
principally clear but rather sophisticated procedure of watching the fate 
of each age cohort that spreads in the set of size cohorts. Let us examine a 
simple example. Consider the fate of a cohort born at time 1, with 
number of trees n(t, 1). At time r+ k there will be n(t + k, k + 1) = gig; 
g,n(r1) trees in size class k+1, and n(t+k,k)= gigs 
Bx-i(1 — gy — dy)n(t, 1) trees in size class k, etc. Summing up all these 
numbers gives us the total number of trees of age t = k: 


a(t+k, n= » n(t + Kk, i). 
t 

It seems to us that base model 2 has some advantages over Eqn 
(13.25). For the latter, in the case of several species, a common set of size 
intervals (Eqn (13.23)) makes the species’ growth resolution different; a 
particular set makes the whole system almost inoperable (e.g. if we need, 
as is often the case, to compare species heights for defining competitive 
relations). This approach is disquieting when one has to deal with a 
number of tree species with noticeably different growth rates and, 
moreover, try to include a description of, sav, grass and shrub growth. 
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System (13.26) does not have any of these disadvantages. The size-class 
approach can, of course, be useful for one-species even-aged tasks aimed 
at analyzing size distribution dynamics. | 

In the next section we will show the simulation abilities of base model 2 
using an example of the modeling of non-stationary age distribution 
behavior. 


13.4.2. Empirical age dynamics in simple forests 


First of all, we will describe the non-stationary age distribution behavior 
of a generalized coniferous species. Such a situation originates after a 
severe catastrophic disturbance, which entirely obliterates the initial 
‘maternal’ stand and provides zero initial conditions 


n(t=0,t) =0 


This situation is typical of major fires in taiga forests, total phytophagy 
defoliation (after which the trees quickly wilt and die). hurricanes, and, 
of course, cutting by man. The subsequent successions have been repeat- 
edly described and analyzed in the literature (Semetchkin 1970; Leak 
1975; Hett & Loucks 1976; Larson & Oliver 1979; Oliver 1981; and 
others). 

Further, the intensity of seed influx f will be the central ecological 
parameter. Clearly, at sufficiently low f values, the area will be inhabited 
by a population with a monotonic age distribution dn/dt « 0. which 
attains equilibrium during one generation time. This case is typical for 
habitats with poor soil and climatic conditions. 

In comparatively better conditions, and when f is large enough, 
pioneer individuals capture the area and exhaust the resources (light, soil 
oxygen, nitrogen, etc.), and the seedlings of the next ages die off owing to 
competitive suppression. The result is a ‘package’ of older trees and a 
zero-gap at that part of the age distribution which corresponds to younger 
trees (Fig. 13.6), a picture that has been described repeatedly (Zubarev 
1965; Semetchkin 1970; Kazimirov 1971; Leak 1975; Sprugel 1976; Oliver 
1978, 1981; Francline & Waring 1979; Larson & Oliver 1979; Glebov & 
Kobyakov 1984). For a full review of different types of non-stationary age 
distributions, see Katayeva & Korzukhin (1987). In all cases shown in 
Fig. 13.6, we have the ‘running wave’ or several waves along the t-axis. 

The subsequent dynamics may go one of two ways: (a) after one or 
several damping waves the age distribution attains equilibrium; or (b) the 
system enters an oscillating regime which can be either fully or quasi- 
periodic. For a full review of field observations, see Katayeva & Korzuk- 
hin (1987). Here we will list only some typical cases. 

A prolonged (100-200 years) endogenous periodic cycle was reported 
for deciduous forests in the far eastern USSR (Kolesnikov 1956; Vasiliev 
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& Kolesnikov 1962; Rozenberg, Manko & Vasiliev 1972), for fir forests in 
the Ural region (Smolonogov 1970; and others), for beech Crimea forests 
(Sukachev & Poplavskaya 1927), and for balsam fir and eastern hemlock 
forests in the Great Lakes region (Hett & Loucks 1976). An endogenous 
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Fig. 13.6." Age waves’ emerging after major forest cover disturbance 
in sufficiently dense stands. n(r, 7) is the number of trees in age 
interval Ar = 1 under age t; t is moment of observation (equal to time 
since disturbance). Cases of strong suppression of juveniles by 
mature canopy: (a) post-fire succession, Pinus sibirica, East Siberia 
(Semetchkin 1970); (b) post-cutting succession. sum of oak, maple, 
and birch, New England (Oliver 1978). Cases of weak suppression: 
(c) post-fire succession, sum of spruce, Pinus sibirica, and fir. West 
Siberia (Korzukhin, Sedych & Ter-Mikaelian 1988); (d) post-fire 
succession. mixed oak forest. Appalachian Mountain (Ross, Sharik 
& Smith 1982). 
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periodic cycle of 60-80 years was reported for Abies forests in 
New England (Sprugel 1976; Reiners & Lang 1979; Foster & Reiners 
1983; Moloney 1986) and for Abies forests in Japan (Oshima er al. 
1958). 

Unfortunately, only a few papers (e.g. Hett & Loucks 1976; Reiners & 
Lang 1979; Kozin 1982) give information about age distributions. Of 
great importance are the data about the size of area which is occupied by a 
single-phase stand: this area is usually between 0.02 and 1 ha (Oshima er 
al. 1958; Foster & Reiners 1983; Moloney 1986). 

A further set of works is dedicated to spatial waves in pure dark- 
needled coniferous stands (Oshima er al. 1958; Sprugel 1976; Reiners & 
Lang 1979; Foster & Reiners 1983; Moloney 1986). If the system is 
observed in a single ‘point’ — reported above to be between 0.02 and 1 ha— 
it exhibits periodicity of age distributions of trees. Spatial waves are 
generated by synchronization of phases in different points owing to 
strong, unidirected winds. When there are no winds the phases of 
development become different. 

The next section is devoted to a modeling analvsis of periodic 
phenomena. 


13.4.3. One-species oscillation models 


First. we offer two simplified qualitative models (for details see Korzuk- 
hin (1980), Antonovski, Korzukhin & Matskiavichus (1989) and Korzuk- 
hin. Antonovski & Matskiavichus (1989)). which are some particular 
cases of the general multi-aged base model in Eqn (13.26). 

The first system is differential and describes number dynamics in a 
system with three age classes (n; are their tree numbers): 


ni = V(yo)f — ni 
ña =n — ns (13.27) 


Az = na — n3- 


Here we ignore the mortality in class 2 and consider the seed background 
as the only source of regeneration: its intensity is f and offspring survival 
is V, depending on total leaf area, g = Sana + San; (Sir, is neglected), 
and y is a generalized interaction parameter (for example, light extinction 
coefficient). By using trivial qualitative equilibrium point analysis, it is 
easy to show that system (13.27) can realize sustainable oscillations. For 
example. when V = exp (—yo), this happens when 


8 
fec 25 P BG 80A, 259]. 
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The discrete analog of system (13.27) 
nlt + 1) — V(ya)f 


n»(t + 1) = ni(t) (13.28) 
nyt 5 1) = ng (t) 
where o(t) = X S(i)n,(t), with stepwise viability function and linear law of 


wt 
ino qot 1 when o < o* 


V(yo) = | _, 5S) Sali- 1) 
0 otherwise 

also realizes sustainable age oscillations. For example. when 

af > o*/(T — 1), there is a regime with period 


BT cSt 


E] 


[z] denotes the integer part of z. The solution of (13.28) with zero initial 
conditions has the form of periodically running waves with K non-zero 
age classes in a package: 
ní(t) = (0.....0)— (f,0.....0)— (f. f, 0... .. 0) = à 
—» (0... S cs oO. 40) (Gass; 0, f) 5 (0.....0). 
The numerical and realistic oscillation model (see the source works 
cited) originates from the general form of base model 2 (Eqn (13.26)). 
Crown morphology was the most simple: crowns were horizontal. and 
again with areas (Eqn (13.10)) 
S, (t. 7) = am^(t. t) (13.29) 
(a and a are parameters); decay and respiration are given by Eqn (13.11). 
Assimilation was proportional to normalized light intensity 
A(Rj) = bRj(t. x) Ryo, A(Rjo) = Amas = b. (13.30) 
where R,is given by (Korzukhin & Ter-Mikaelian (1982); cf. Eqn (13.13)): 


R(t, t) = Ryo sin v exp [uas T)n(t, T) 


=% Fi SL int) 
u-rl 
q is the angle of the sun. Within total viability, Y = WV, the density- 
independent component was constant: 


W(t) = Wo «€ 1, (13.31) 
and the density-dependent one had the form 
1— vf A(R;) 
1—4 Iz] h jl UI. 
V(v) = pea] UM dw (tS 


0 otherwise 
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(q 5 1, u< 1, and B are parameters), that realizes the possibility to put 
smooth (small £) and near-stepwise (large £) viability functions (Fig. 
13.1). Only the background seed quantity was considered: as in Eqn 
(13.28), p(t) = 0. 

Finally, we dealt with the system 


n(t+1,1)=f 
n(t+1,7+1) 2 WoV(v,) (13.33) 
m(t + 1,1) = mo z1.::41 


m(t+ l1, r 4 1) 2 m(t, 1) + S; (m,)A(R,) — cm(t, v) 
[m, = m(t, t), R, = R(t, T), v, = A(R,)/A max] 
where V, S, A, and R are the functions defined above. 


Parameter definition 


Species lifespan T was set to 200 years. The value of a in Eqn (13.29) most 
commonly lies in the interval 0.5 = a = 0.8; we took a = 0.7. At the end 
of the lifespan of coniferous species, the ratio a, = m7" /mg,, = 0.05 
(Larcher 1975), where m, is leaf biomass. The mass-surface transition 
factor for needles was taken as a, = m, /s; = 0.15 [kg fresh mass wh 
After taking Mmax = 1000, we can calculate that a = aym\,%/a, = 2.65; 
this means S}** = 334 m? and mP™ = 50 kg. Ry can be adopted as 1; 9 
equaled 0.76 [rad). A typical value of b in Eqn (13.30) is in the order of 
10^! [kg fresh mass m ^ yr^!]; we took b = 0.15. The parameter c of 
unified respiration and decay losses can be found from the condition 
m(t)— Myx when t—°% and competition is absent: this gives 
amS axb sin @ — CMmax = 0 or c = 0.034 [kg ' yr]. W in Eqn (13.31) was 
equal to 0.98053 (by the age T = 200, 2% of trees are left); q equaled 0.5, 
that is, a moderate jump in V at z — u was allowed. 

The values of y, f, B and u were varied in the search for oscillations. In 
this way, we experimented with generalized coniferous tree species 
having plausible parameter values from measurements and a number of 
free parameters for the searching of oscillations. Among the latter 
y = 0.3-0.8, f usually varies from some hundreds to thousands, u is near 
the light compensation point for the whole tree, that is u = 0.1-0.3, and B 
measures the plasticity of tree response to shortage of resources, a value 
which is unknown. The initial conditions comprise zero densities, 
n(0, rt) 2 0 and mo = m(0, 1) = 107? [kg]; m(0, t) = 0 for cz 2. After 
many numerical experiments we have concentrated on the value f — 60, 
which relates to all the adduced results. The system calculated up to a /7** 
of 2000 years. 

All dynamical regimes in three-dimensional parametric space (y, f, u) 
are given in Table 13.1, where the period of oscillations, 0, when they 


362 M. D. Korzukhin and M. Ya. Antonovski 


Table 13.1. Oscillation period 0 (y. f. u) and the type of dynamical 
behavior in Eqn (13.33) 


RO, 86 
RO, 82 
QO. 52-68 


RO, 86 
QO, 43-71 
QO, 49-66 


Note: SE = stable, RO = regular, and QO = quasi-regular oscillations. 


exist, is also shown. When the oscillations are quasi-regular. we adduce 
the observed interval of 0. It can be seen that low values of y, f, and u 
promote a stable behavior: as their values increase the system becomes 
unstable and finallv falls into a stochastic regime. 

The behavior of some important characteristics as functions of y for 
parameter values f= 1000. u — 0.15 is presented in Fig. 13.7. Up to 
y = 0.2, the system attains stable equilibrium. Near y = 0.22 it becomes 
unstable. and a stable periodic regime arises (Fig. 13.8). Figure 13.7 
shows the dependence on y of oscillation amplitude n™* (r, 100) and n™" 


(r. 777), where 7" is the minimum age when the maximum number is at 


7™* = 100 (see also Fig. 13.8 where curves n(t. t) are pictured at the time 
when the maximum goes through age 100). Note (Fig. 13.7) the high 
stability of the oscillation period, 80 = 8 = 88. and the value of n™** 


whereas n™" is greatly reduced. 


Leaf area index 

The time behavior of leaf area index, o(t) = X,S(t. tT)n(t. t), shows the 
appearance of expected dynamical invariant (13.8): bevond r = 50-100, a 
assumes a relatively constant value (the relative variations 
Agíg^ *'*8* = (.15) which is the maximum possible value for the given site 
and species parameters. This behavior corresponds well with the situ- 
ations observed in even- and multi-aged stands. At y = 0 (Fig. 13.7) (free- 
growing trees) o™* is very large and equals 115, a totally unrealistic value 
for natural populations. However, o(y) then falls quickly, attaining 
reasonable values by the time y is only 0.1. 


Ecophysiological invariant 

Among the effects discovered the greatest interest was aroused by the 
appearance of some parametric invariant. In our previous works (Anto- 
novski, Korzukhin & Matskiavichus 1989; Korzukhin. Antonovski & 
Matskiavichus 1989) we have found it occasionally. but now we are able 
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Parameter magnitude 


0.1 02 0.4 0.6 0.8 1.0 


q 
Fig. 13.7. Some quantities of system (13.33) with various values of 
the parameters and extinction coefficient y. Two curves (filled 
circles) arising after y 0.22 show amplitude of the younger wave; for 
meanings of n™* = n(t, 100) and n™ = n(t,7™") see Fig. 13.8. Open 
circles, oscillation period 8 (y); squares, leaf area index 6 (y) in T= 
2000; pluses, value of J (13.35). 


to develop some theoretical base. Let us put the following question: what 
is the interrelation between maximum leaf area index, o0™**, which a 
stand develops, extinction coefficient, y, and leaf (or tree) viability? It is 
clear that when y is large, under fixed leaf shade tolerance a stand can 
sustain a fewer number of layers than when y is small. Physiologists know 
well that the lowest leaf layer in a stand exists only when its respiration 
losses are less than net assimilation. Below, we try to express this rule in 
model terms. 
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Let us consider the ideal situation when viability function V is step- 
wise: 
1, when Rj/Rjy > u 
WR) = . 
0 otherwise; 


that is, a tree with height A dies when for relative light supply the 
following condition is true 


R,(h)/Ryo S u. (13.34) 
Further, when a stand develops ø leaf layers, the last layer receives 
exp (—yo) 


share of initial light. So, the condition of existence of exactly d™™* layers, 
obtained from Eqn (13.34), is 
exp (—yo""*) = u. 

It is suitable to introduce the function 

J(o. y.u) = u exp (y0). ] (13.35) 
which can be named as an 'ecophysiological invariant’. It connects three 
fundamental magnitudes relating to stand functioning. Among them, u 
can be identified with a compensation point on the leaf light curve 


measured for 24-hour leaf carbon balance (normalized to maximum 
light). 


Woe || 


1000 


7 WANS 
: | NE LT 


Fig. 13.8. Age distributions of tree numbers n(t,r): solution of 
system (13.33) under various y. 1, stable distribution for y = 0, n(7) = 
1000 exp( —0.01977); 2,3, 4, stable oscillations for y = 0.226, 0.4, 0.8; 
curves are taken in time moment t! when maximums were going 
through age — 100; see text for parameter values. 
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Table 13.2. Dependence of leaf area index o (first number) and 
invariant J (second number) on system parameters ` 


For our simplified viability model the equality 

J(o, y, u) =J™™* =1 
should be valid. In reality, many factors deviate J from 1. Any tree has a 
finite probability of dying even if its relative light supply is greater than u. 
So, a stand cannot reach the value o™*, which corresponds to the 
stepwise form of V; and one will expect g < o™* and the value of J will be 
less than 1. Deviation of J from 1 is a measure of stand structure deviation 
from optimum structure. The concrete value of J is defined by the form of 
the V function and other system parameters. Returning to our model, we 
adduce the example of J calculation at u = 0.15 (Fig. 13.7) and for u = 0.1 
and u = 0.2 (Table 13.2). Certain scattering and non-monotonic behavior 
of J can be explained by taking its values at ¢ = 2000 and in an oscillation 
regime. 

It seems to us that further development of this approach lies in 
calculations of ecophysiological invariant (13.35) for real stands, while of 
central interest is watching the degree of J stability with parameter 
variations and time, and analysis of its deviations from value J™* = 1. 

Table 13.2 also shows an approximate constancy of g as a function of 
seed migration, f, under a given y and u. Here we are witnessing a 
‘habitat-saturation’ effect, which is observed in real ecosystems and 
explained, as above, by maximum light resource utilization. The effect 
was also studied by Korzukhin (1986). 

Numerical experiments showed that a large curvature of viability 
function V(v) Eqn (13.32) is essential for the appearance of the oscil- 
lations. The influence of q, u, and f can be seen from the value of the 
derivative at v = u: V’(v)|,-, = —qB/(1 — u). It is worth mentioning that 
the used value of B = 60 yields a ‘near-step’ function V(v), which is close 
to purely stepwise viability functions widely used in forest models (see 
review in $13.2.2). 

The dynamics of formation and passing of the first age wave is of great 
interest for forest ecology regardless of the question of oscillations. The 
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n(t, 7) 
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Fig. 13.9. Solution of system (13.33) under y = 0.4 and u ="0.2: 1. 

first ‘age package’; 2, stable oscillations. 
first wave is much higher and narrower than the succeeding ones (Fig. 
13.9). One can see that at z = 100, the system is on the edge of the next 
wave origin. 


13.5. From stand to forest; addition of above-ground species 


Up to this stage, we have discussed stand dynamics only, which is a 
serious contraction of real forest composition. However, although trees 
can dominate in certain senses, e.g. in terms of the proportion of live 
ecosystem biomass, essentially they can depend on other species. Among 
these, the above-ground plant species such as grasses, mosses and lichens 
are most important for boreal forests because they provide the boundary 
conditions for tree regeneration. In this aspect, trees as a life form exhibit 
an apparent weakness in their regeneration strategy: many trees suppress 
mosses and grasses and simultaneously their own seedlings. Sparse stands 
allow seedlings to grow but also promote their competitors. The interplay 
between these extremes can lead to interesting dynamical and life- 
strategy scenarios. Below, we describe some of the model embodiments 
of a simplified ecosystem consisting of one tree and one above-ground 
species as an example of a two-life-form community. Modeling activity in 
this field is very limited compared with tree systems only, so we can 
adduce only a few results. 


13.5.1. Modeling of moss dynamics 


Many dark-needled coniferous boreal forests are characterized by a 
noticeable moss organic layer on the forest floor. The thickness of this 
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laver, taking into account both live and dead parts, is up to between 30 
and 50 cm. The layer is an important structural component of a forest, 
controlling energy flow, nutrient cycling, water relations, and, through 
these, stand productivity and dynamics (Bonan & Shugart 1989). For 
example, soil temperature and depth of permafrost are directly related to 
the thickness of the layer. Another example is the dependence of the 
viability of different types of tree seeds on moss layer thickness (Fig. 
13.10), which will be discussed below. 

First, let us consider only live moss. The simplest formalization of its 
growth can be made in terms of carbon balance using one variable: green 
moss biomass u [kg m7]: 

A m S,A,C(R,) — €H. (13.36) 
where respiration and decay are assumed to be proportional to uw, R, 
denotes mean light falling per unit area. A,, is specific assimilation. and 
leaf area S, can be taken as proportional to 4 (S, = a4). The central idea 
of this model is that moss growth must be auto-restricted, that is. Eqn 
(13.36) must give 4 — u ™®™ in the absence of trees and other competitors. 
If we consider light extinction in a relatively thin moss layer with 
thickness / and a vertical profile of light R,,(/), then. approximately. 


R, = R(h/2) = Ry exp (775,4/2). 


where A, is the light at the top of the layer. If we suppose that A, is 
proportional to R, (A,, = bR,,), this gives us a simple equation of moss 
layer growth: 


f= a,uUbR, o exp (—ya,u/2) — cu = v exp (—qu) — cu. (13.37) 


Ulu) 


u 


Fig. 13.10. Qualitative dependence of seed survival and germi- 
nation, U, on moss biomass (or thickness) 4 for boreal tree species. 
(1) Small seeds (birch, spruce); (2) large seeds, in particular Pinus 
sibirica, whose seed viability is maximal under a certain moss 
thickness (summarized in Katayeva & Korzukhin (1987)). 
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For sustainable existence of moss, v must be greater than c. This model 
coincides with that proposed for grasses by McMurtrie & Wolf (1983). A 
more realistic model was developed by Bonan & Korzukhin (1989) where 
the second moss variable v (dead moss biomass [kg m ?]) and also the 
usual light curve for assimilation were taken into account: 

jn 8, ee cp 
v Ct — C3V. 


In this equation, c, is the specific decay rate of live moss. c» is its 
respiration losses, c, is the specific oxidation rate of dead moss. a, and a> 
are light curve parameters, and a; is the compensation point. The model 
was identified with the help of actual ecophysiological (Larcher 1975) and 
field (Van Cleve et al. 19835) moss data. 


13.5.2. Lichen-trees system modeling 


The principal components for the construction of various models of the 
discussed type have been identified in previous parts of this chapter. For 
example. we may join the age-distributed base tree model 2 (Eqn (13.26)) 
with the simplest moss growth equation (Eqn (13.37)). In fact. this 
equation can be applied to any above-ground species because principles 
on which the equation is based are very general. For the sake of diversitv, 
let us consider the lichen-trees situation that is described in detail by 
Sirois in Chapter 7 of this book. Among other properties of the system. 
there is a possibility of sustainable occupation of a given territory either 
by trees (represented by spruce) or by lichen: this depends on initial 
conditions of area settlement. 

In the simplest case, lichen-tree interactions are expressed (a) in terms 
of dependence of incident light reaching lichen layer. R;;, on tree leaf 
area a: Ro = R;o(0), and (b) in terms of seedling viability dependence on 
lichen thickness, U = U(4). These arguments, 4 and g, interconnect tree 
and lichen dynamics. 

Whereas the complete model (Eqns (13.26) and (13.37)) is difficult to 
analyze, some qualitative properties can be demonstrated on the basis of 
its simplified version. Let us take one generalized variable N (total 
number of trees) and consider only seed background. A very simple 
lichen-tree model turns out to be 

B MOM R (13.39) 
A = [v exp (CkN — qå) —c]A. 
where K is proportional to average single-tree leaf area multiplied by the 
extinction coefficient for trees. In this approximate description, N can be 
changed to total tree biomass M. The behavior of this system is typical of 
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two-dimensional systems with ‘competitive’ interactions, which are 
encountered repeatedly in forest modeling (e.g. McMurtrie & Wolf 
1983). 

First consider the case of the suppression form of function U(A) only 
(Fig. 13.10). Let us take U(4) 2 1 — A/4* when 4 € 4*; otherwise 
U(A) = 0 (A* is a border lichen thickness). Simple qualitative analysis of 
system (13.39) gives us four standard situations. By denoting r = In (v/c), 
we have: 

1. If kf ^ rand q4* > r, we have one stable equilibrium (f, 0), that is. 
trees outcompete lichen. 

2. If kf « r and qà* >r, we have one stable equilibrium (N°, 2°), 
where N° < f, 2° < 4*, that is, lichen and trees coexist; N° and 2" can be 
found from the linear system that originates from Eqn (13.39) under N, 
À — 0: 


kN -qi =r; N — (1 — A/4*)f = 0. 


3. If Af>r and q4* «€r, we have two stable equilibria (f, 0) 
and (0. 7/g). which are realized depending on initial conditions N(0). 
A(0). 

4. If kf < rand qå” < r, we have one stable equilibrium (0. r/q). that 
is, lichen outcompetes trees. 

The case of a non-monotonic function U(Z) (Fig. 13.10) adds the 
possibility of the coexistence of two stable states: (Nj, 2°) and (NÌ. 0). 


Now let us imagine a movement along an environmental gradient from 
dry and warm conditions, which are favorable for tree growth and 
unfavorable for lichen (fis large, v is small), to wet and cold conditions 
where the reverse is true. In this case, system (13.39) will realize the 
consequences of the stable state moving from 1— 2— 4 or 1— 3— 4, 
that is, the transition from temperate forest without lichen to tundra 
without trees, passing between these extremes through the possibility of 
lichen-tree coexistence (boreal forest) (Fig. 13.11). There is absolutely 
no doubt that any generalization of this qualitative system (e.g. Eqns 
(13.26) and (13.38)) will retain this transition effect. 

In a more realistic imitative, but not mathematical, model (Bonan & 
Korzukhin 1989) a very similar type of behavior was actually discovered 
in the field for tree-moss communities disposed along a cold—wet to hot- 
dry environmental gradient. The sites were typical of central Alaska and 
were dominated by black and white spruce, white birch and trembling 
aspen, plus sphagnum mosses. The two-compartment moss model (Eqn 
(13.38)) was used. Functions of seed germination suppression, U,(M) for 
each tree species i, were reconstructed from field data; here M, = u + vis 
the total moss biomass. The dynamics of the trees were described by using 
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(a) 


DAWs) 
f(Ws) 
(b) Ws 
Temperate | Boreal | Tundra 
forest forest I 
| 
| 
I 
N9 (Ws) l 
M (Ws) 


N° [p(W 5), f (W 9] A9 [p(W s), FWs)] 


state 1 state 2 state 4 Ws 


Fig. 13.11. Parametric evolution in lichen-tree system (13.39): if one 
considers p and f as monotonic functions of soil moisture Wg (a). then 
steady states of (13.39). N"(W.). 2"(W;). realize trajectories (b) 
which can be interpreted as different types of lichen-tree ecosystems 
disposed along an environmental gradient (a case of monotonic L(A). 
as a consequence of states: | — 2 — 4). The conditions of realization 
of these states in system (13.39) are as follows. 

State 1: In [p(WJ/c)] < kf(Ws) <q* , 

State 2: Af(Ws) < In [Cp(Ws/c)] < q^ . 

State 4: kf(Ws) < q^* < In [p( W y/c)]. 


a common gap-modeling technique. Field observations along the gradi- 
ent supported model results; see the cited work for details. 


6. Concluding remarks 


In this chapter, two aspects of forest modeling were touched upon. The 
first relates to the whole field of stand modeling and reviews the 
compulsory, elementary blocks which form the base of every type of 
model. Here, we presented the most contemporary methods of tree 
growth, inter-tree competition, and tree viability formalization. 

The second modeling aspect concerned developing and running 
models of forest dynamics at the population level. These models are 
based on the use of certain dynamical equations derived with the help of 


mechanistic considerations. We demonstrated the researching abilities of 
these models by means of explanation of some central dynamical effects 
observed in even- and multi-aged forest stands. Among these effects, 
wave behavior of tree age distributions was the most demonstrative 
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example of the abilities of population-level equations. 
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preparing the text for publication, and, of course, for her benevolent patience 
during numerous modifications of the paper. The authors would also like to 
express thankfulness to their colleagues in forest modeling for beneficial dis- 


cussions. 


Appendix 


List of variables and parameters used 


specific assimilation rate [kg m? yr—'] 
tree mortality (21 — Y) [yr !] 

tree diameter at breast height [cm] 

tree decay losses [kg yr !] 

seed immigration rate [ha ^! yr" !] 

stand cross-sectional area (=D*N/4) [m°] 
light extinction coefficient [dimensionless] 
tree height [m] 

lichen biomass [kg m7] 

total tree biomass [kg of dry or fresh mass] 
tree leaf biomass [kg] 

total moss biomass (=u + v) [kg m7] 
green moss biomass [kg m^ ?] 

tree population age distribution [ha 
tree population density [ha !] 
ecophysiological invariant [dimensionless] 

dead moss biomass [kg m] 

tree seed production [yr !] 

population seed production [yr^ !] 

light intensity [arbitrary units] 

incident light on moss and lichen surface [arbitrary units] 
vertical profile of light in moss layer [arbitrary units] 
amount of soil resource [arbitrary units] 

tree respiration losses [kg yr *] 

tree leaf, root, and crown projection area [m7] 

moss leaf area [m7] 

population leaf area index [dimensionless] 

population crown projection area [dimensionless] 


yr] 
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population leaf area index above height A [dimensionless] 
time [yr] 

tree age [yr] 

seed viability connected with moss and lichen effect [yr !] 
density-dependent component of tree survival probability 
[yr] 

density-independent component of tree survival 
probability [yr !] 

tree survival probability (—1 — d) [yr^ !] 


